/*	*************************************************************/
/*     	File Name:	Final do Models for RR.do				   */
/*     	Date:   	December 1, 2023                          */
/*      Author: 	Kellan Borror                              */	
/*      Input Files: Final Data.csv	    						*/
/*     	Output File: 						                    */	
/*	*************************************************************/



**Set Working Directory and import data**
cd "/Users/kellanrobinson/Desktop/Final Paper and Analysis for Submission 1"
import delimited "/Users/kellanrobinson/Desktop/Final Paper and Analysis for Submission 1/Final Data.csv"

**generate logged variables**
gen log_force_size = ln(force_size)
gen log_pop_density = ln(pop_density)
gen log_gni_per_cap = ln(gni_per_cap)
gen log_gdp_mission = ln(gdp_mission)
gen log_pko_density = ln(force_density)


**Models with Zero-inflated Binomial Regression**
**DV= sea_cat1_total, IV= log_ngo_counts**

**Model 1**
zinb sea_cat1_total log_ngo_counts pts_mission gender_ratio log_force_size log_pop_density log_gni_per_cap, inflate (log_ngo_counts pts_mission gender_ratio log_force_size log_pop_density log_gni_per_cap) cluster(mission)

**Model 2**
zinb sea_cat1_total log_ngo_counts corruptionmilvars gender_ratio log_force_size log_pop_density log_gni_per_cap pko_fatalities_malicious, inflate (log_ngo_counts corruptionmilvars gender_ratio log_force_size log_pop_density log_gni_per_cap pko_fatalities_malicious) cluster(mission)

**Model 3**

zinb sea_cat1_total log_ngo_counts log_gdp_mission gender_ratio log_force_size log_pop_density log_gni_per_cap, inflate (log_ngo_counts log_gdp_mission gender_ratio log_force_size log_pop_density log_gni_per_cap) difficult cluster(mission)

**Model 4**

zinb sea_cat1_total log_ngo_counts pts_mission gender_ratio log_pop_density force_density, inflate (log_ngo_counts pts_mission gender_ratio log_pop_density force_density) difficult cluster(mission)

**Model 5** 

zinb sea_cat1_total log_ngo_counts corruptionmilvars force_qual_mil gender_ratio log_pop_density force_density, inflate (log_ngo_counts corruptionmilvars force_qual_mil gender_ratio log_pop_density force_density) difficult cluster(mission)

**Model 6**

zinb sea_cat1_total log_ngo_counts pts_mission force_qual_mil gender_ratio log_force_size ln_km2_country , inflate (log_ngo_counts pts_mission force_qual_mil gender_ratio log_force_size ln_km2_country) cluster(mission)


**Model 7**

zinb sea_cat1_total log_ngo_counts corruptionmilvars gender_ratio log_force_size ln_km2_country, inflate (log_ngo_counts corruptionmilvars gender_ratio log_force_size ln_km2_country) difficult cluster(mission)


**Models with Zero-inflated Binomial Regression**
**DV= sea_cat1_total, IV= log_ngo_counts_hr**

**Model 8**
zinb sea_cat1_total log_ngo_counts_hr pts_mission gender_ratio log_force_size log_gni_per_cap, inflate (log_ngo_counts_hr pts_mission gender_ratio log_force_size log_gni_per_cap) cluster(mission)

**Model 9**
zinb sea_cat1_total log_ngo_counts_hr corruptionmilvars gender_ratio log_force_size log_pop_density log_gni_per_cap pko_fatalities_malicious, inflate (log_ngo_counts_hr corruptionmilvars gender_ratio log_force_size log_pop_density log_gni_per_cap pko_fatalities_malicious) cluster(mission)

**Model 10**

zinb sea_cat1_total log_ngo_counts_hr log_gdp_mission gender_ratio log_force_size log_pop_density ln_km2_country, inflate (log_ngo_counts_hr log_gdp_mission gender_ratio log_force_size log_pop_density ln_km2_country) difficult cluster(mission)

**Model 11**

zinb sea_cat1_total log_ngo_counts_hr pts_mission force_qual_mil gender_ratio log_force_size log_gni_per_cap, inflate (log_ngo_counts_hr pts_mission force_qual_mil gender_ratio log_force_size log_gni_per_cap ) difficult cluster(mission)



**Models with Zero-inflated Binomial Regression**
**DV= sea_total, IV= log_ngo_counts_w**

**Model 12**
zinb sea_total log_ngo_counts_w pts_mission  gender_ratio log_force_size log_pop_density log_gni_per_cap , inflate (log_ngo_counts_w pts_mission gender_ratio log_force_size log_pop_density log_gni_per_cap  ) difficult cluster(mission)

**Model 13**
zinb sea_total log_ngo_counts_w corruptionmilvars gender_ratio log_force_size log_pop_density pko_fatalities_malicious log_gni_per_cap, inflate (log_ngo_counts_w corruptionmilvars gender_ratio log_force_size log_pop_density pko_fatalities_malicious log_gni_per_cap) difficult cluster(mission)

**Model 14**

zinb sea_total log_ngo_counts_w log_gdp_mission gender_ratio log_force_size log_pop_density , inflate (log_ngo_counts_w log_gdp_mission gender_ratio log_force_size log_pop_density ) difficult cluster(mission)

**Model 15**

zinb sea_total log_ngo_counts_w pts_mission force_qual_mil gender_ratio log_force_size ln_km2_country log_gni_per_cap , inflate (log_ngo_counts_w pts_mission force_qual_mil gender_ratio log_force_size ln_km2_country log_gni_per_cap ) cluster(mission)





**Apendix** 

**Mil Models** 

**Model 16**

zinb cat1_sea_mil log_ngo_counts pts_mil gender_mil log_force_size log_pop_density log_gni_per_cap, inflate (log_ngo_counts pts_mil gender_mil log_force_size log_pop_density log_gni_per_cap) cluster(mission)


**Model 17**

zinb cat1_sea_mil log_ngo_counts_hr pts_mil gender_mil log_force_size log_gni_per_cap, inflate (log_ngo_counts_hr pts_mil gender_mil log_force_size log_gni_per_cap) difficult cluster(mission)

**Model 18**

zinb sea_mil log_ngo_counts_w pts_mil  gender_mil log_force_size, inflate (log_ngo_counts_w pts_mil gender_mil log_force_size ) difficult cluster(mission)

**pol Models**

**Model 19**

zinb cat1_sea_pol log_ngo_counts pts_mission_pol gender_pol log_force_size log_pop_density , inflate (log_ngo_counts pts_mission_pol gender_pol log_force_size log_pop_density) difficult cluster(mission)

**Model 20**

zinb cat1_sea_pol log_ngo_counts_hr pts_mission_pol  log_force_size, inflate (log_ngo_counts_hr pts_mission_pol  log_force_size) difficult cluster(mission)

**Model 21**

zinb sea_pol log_ngo_counts_w pts_mission_pol   log_force_size log_pop_density, inflate (log_ngo_counts_w pts_mission_pol  log_force_size log_pop_density) difficult cluster(mission)

********************
*** Figures ***
********************

********************
*** Figure 3 ***
********************

*** Run model and save estimates ***
eststo m1: zinb sea_cat1_total log_ngo_counts pts_mission gender_ratio log_force_size log_pop_density log_gni_per_cap, inflate (log_ngo_counts pts_mission gender_ratio log_force_size log_pop_density log_gni_per_cap) cluster(mission)

*** Loop to set global macros that hold natural logged counts of NGOs ***
* Since range is 7 - 330, I will do 10 - 330
macro drop _all
forvalues i = 1(1)33 {
	global set_`i' = ln(`i'0)
}


*** Estimate predicted counts based on set values above ***
estimates restore m1
margins, at(log_ngo_counts = ($set_1 $set_2 $set_3 $set_4 $set_5 $set_6 $set_7 $set_8 $set_9 $set_10 $set_11 $set_12 $set_13 $set_14 $set_15 $set_16 $set_17 $set_18 $set_19 $set_20 $set_21 $set_22 $set_23 $set_24 $set_25 $set_26 $set_27 $set_28 $set_29 $set_30 $set_31 $set_32 $set_33 )) atmeans post saving(mars, replace)


*** Import predicted values *** 
use mars.dta, clear


*** Keep point estimates, upper bound, and lower bound ***
keep _margin _ci_lb _ci_ub


*** Generate counts of NGOs for the graph ***
gen ngo_counts = [_n] * 10


*** Make a graph ***
graph twoway line _margin ngo_count, clwidth(medium) clcolor(black) clpattern(solid) sort || line _ci_lb ngo_count, clpattern(dash) clwidth(thin) clcolor(black) sort || line _ci_ub ngo_count, clpattern(dash) clwidth(thin) clcolor(black) sort legend(off) xlabel(10(40)330) xtitle("Number of NGOs", size(medsmall)) ytitle("Predicted Count of Category 1 and SEA Misconduct", size(medsmall)) scheme(sj) title("Effect of NGO Counts on Category 1 and SEA Misconduct", size(medlarge))
graph export "Model_1.jpg", as(jpg) name("Graph") quality(90) replace


********************
*** Figure 4 ***
********************

*** Import pkat .csv file ***
import delimited "Final Data.csv", clear  


*** Generate logged variables ***
gen log_force_size = ln(force_size)
gen log_pop_density = ln(pop_density)
gen log_gni_per_cap = ln(gni_per_cap)
gen log_gdp_mission = ln(gdp_mission)


*** Run model and save estimates ***
eststo m1: zinb sea_cat1_total log_ngo_counts_hr pts_mission gender_ratio log_force_size log_gni_per_cap, inflate (log_ngo_counts_hr pts_mission gender_ratio log_force_size log_gni_per_cap) cluster(mission)

*** Loop to set global macros that hold natural logged counts of NGOs ***
* Since range is 9 - 97, I'll do 10 - 100 
macro drop _all
forvalues i = 10(5)100 {
	global set_`i' = ln(`i')
}

*** Estimate predicted counts based on set values above ***
estimates restore m1
margins, at(log_ngo_counts = ($set_10 $set_15 $set_20 $set_25 $set_30 $set_35 $set_40 $set_45 $set_50 $set_55 $set_60 $set_65 $set_70 $set_75 $set_80 $set_85 $set_90 $set_95 $set_100 )) atmeans post saving(mars2, replace)


*** Import predicted values *** 
use mars2.dta, clear


*** Keep point estimates, upper bound, and lower bound ***
keep _margin _ci_lb _ci_ub


*** Generate counts of NGOs for the graph ***
gen ngo_counts = 10
replace ngo_counts = 15 in 2
replace ngo_counts = 20 in 3
replace ngo_counts = 25 in 4
replace ngo_counts = 30 in 5
replace ngo_counts = 35 in 6
replace ngo_counts = 40 in 7
replace ngo_counts = 45 in 8
replace ngo_counts = 50 in 9
replace ngo_counts = 55 in 10
replace ngo_counts = 60 in 11
replace ngo_counts = 65 in 12
replace ngo_counts = 70 in 13
replace ngo_counts = 75 in 14
replace ngo_counts = 80 in 15
replace ngo_counts = 85 in 16
replace ngo_counts = 90 in 17
replace ngo_counts = 95 in 18
replace ngo_counts = 100 in 19


*** Make a graph ***
graph twoway line _margin ngo_counts, clwidth(medium) clcolor(black) clpattern(solid) sort || line _ci_lb ngo_counts, clpattern(dash) clwidth(thin) clcolor(black) sort || line _ci_ub ngo_counts, clpattern(dash) clwidth(thin) clcolor(black) sort legend(off) xlabel(10(10)100) xtitle("Number of Human Rights NGOs", size(medsmall)) ytitle("Predicted Count of Category 1 and SEA Misconduct", size(medsmall)) scheme(sj) title("Effect of Human Rights NGO Counts on Category 1 and SEA Misconduct", size(medium))
graph export "Model_7.jpg", as(jpg) name("Graph") quality(90) replace


*********************
*** Figure 5 ***
*********************

*** Import pkat .csv file ***
import delimited "Final Data.csv", clear  


*** Generate logged variables ***
gen log_force_size = ln(force_size)
gen log_pop_density = ln(pop_density)
gen log_gni_per_cap = ln(gni_per_cap)
gen log_gdp_mission = ln(gdp_mission)


*** Run model and save estimates ***
*
* 
*
eststo m1: zinb sea_total log_ngo_counts_w pts_mission  gender_ratio log_force_size log_pop_density log_gni_per_cap , inflate (log_ngo_counts_w pts_mission gender_ratio log_force_size log_pop_density log_gni_per_cap  ) difficult cluster(mission)


*** Loop to set global macros that hold natural logged counts of NGOs ***
* Since range is 1 - 49, I'll do 1 - 50
macro drop _all
forvalues i = 5(5)50 {
	global set_`i' = ln(`i')
}
global set_1 = ln(1)


*** Estimate predicted counts based on set values above ***
estimates restore m1
margins, at(log_ngo_counts = ($set_1 $set_5 $set_10 $set_15 $set_20 $set_25 $set_30 $set_35 $set_40 $set_45 $set_50 )) atmeans post saving(mars3, replace)


*** Import predicted values *** 
use mars3.dta, clear


*** Keep point estimates, upper bound, and lower bound ***
keep _margin _ci_lb _ci_ub


*** Generate counts of NGOs for the graph ***
gen ngo_counts = 0
replace ngo_counts = 5 in 2
replace ngo_counts = 10 in 3
replace ngo_counts = 15 in 4
replace ngo_counts = 20 in 5
replace ngo_counts = 25 in 6
replace ngo_counts = 30 in 7
replace ngo_counts = 35 in 8
replace ngo_counts = 40 in 9
replace ngo_counts = 45 in 10
replace ngo_counts = 50 in 11


*** Make a graph ***
*
* NOTE: I called the first tick a 0 even though you can't use the natural log on a 1. It's for graph continuity. Feel free to change it 
*
graph twoway line _margin ngo_counts, clwidth(medium) clcolor(black) clpattern(solid) sort || line _ci_lb ngo_counts, clpattern(dash) clwidth(thin) clcolor(black) sort || line _ci_ub ngo_counts, clpattern(dash) clwidth(thin) clcolor(black) sort legend(off) xlabel(0 5 10 15 20 25 30 35 40 45 50) xtitle("Number of Women's Rights NGOs", size(medsmall)) ytitle("Predicted Count of SEA Misconduct", size(medsmall)) scheme(sj) title("Effect of Women's Rights NGO Counts on Category 1 and SEA Misconduct", size(medium))
graph export "Model_10.jpg", as(jpg) name("Graph") quality(90) replace




